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Abstract. The effects of relativistic dynamics and thermodynamics in the development of Kelvin-Helmholtz in- 
stabilities in planar, relativistic jets along the early phases (namely linear and saturation phases) of evolution has 
been studied by a combination of linear stability analysis and high-resolution numerical simulations for the most 
unstable first reflection modes in the temporal approach. Three different values of the jet Lorentz factor (5, 10 
and 20) and a few different values of specific internal energy of the jet matter (from 0.08 to 60. Oc^) have been 
considered. Figures illustrating the evolution of the perturbations are also shown. 

Our simulations reproduce the linear regime of evolution of the excited eigenmodes of the different models with 
a high accuracy. In all the cases the longitudinal velocity perturbation is the first quantity that departs from 
the linear growth when it reaches a value close to the speed of light in the jet reference frame. The saturation 
phase extends from the end of the linear phase up to the saturation of the transversal velocity perturbation 
(at approximately 0.5c in the jet reference frame). The saturation times for the different numerical models are 
explained from elementary considerations, i.e. from properties of linear modes provided by the linear stability 
analysis and from the limitation of the transversal perturbation velocity. The limitation of the components of the 
velocity perturbation at the end of the linear and saturation phases allows us to conclude that the relativistic 
nature of the flow appears to be responsible for the departure of the system from linear evolution. 
The high accuracy of our simulations in describing the early stages of evolution of the KH instability (as derived 
from the agreement between the computed and expected linear growth rates and the consistency of the saturation 
times) establishes a solid basis to study the fully nonlinear regime, to be done elsewhere. The present paper also 
sets the theoretical and numerical background for these further studies. 
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1. Introduction 

Relativistic jets are found in many astrophysical ob- 
jects like active galactic nuclei (see eg. Zensus 1997), mi- 
croquasars (Mirabel & Rodriguez 1999) and gamma-ray 
bursts (Sari, Piran & Halpcrn 1999). The basic physical 
property of jets is their huge flux of kinetic energy, taken 
away from compact unresolved or marginally resolved cen- 
tral regions. The presence of jets in such objects provides 
a unique opportunity to study the physics of these central 
regions through investigations of the jet flows. However, 
the extraction of information about the central object 
is not a simple task, because jets interact in complex 
ways with the surrounding interstellar and intergalactic 
medium. On one hand the dynamical interaction of the 
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jet matter with the ambient medium leads to the forma- 
tion of shocks, turbulence, acceleration of charged par- 
ticles and subsequent emission of a broad-range electro- 
magnetic radiation, which makes jets observable. On the 
other hand the complex nature of the interaction which 
is due to the presence of flow instabilities, especially the 
Kelvin-Helmholtz (KH) instability, makes it difflcult to 
distinguish those properties which are directly related to 
the central object from those which are due to the inter- 
action of jets with the ambient medium. As an example 
one can invoke the wavelike patterns in jets, which may 
result from the precession of the rotation axis of the ac- 
cretion disk at the jet base, or from KH instability, which 
finds favorable conditions at the interface of jet and ex- 
ternal medium (Trussoni, Ferrari & Zaninetti, 1983). Very 
recently, the KH instability theory has been successfully 
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used to interpret the structure of the pc-scale jet in the 
radio source 3C273 (Lobanov & Zensus 2001). 

At kiloparscc scales, the surprisingly stable propaga- 
tion of relativistic jets in some sources (e.g., Cyg A) con- 
trasts with the deceleration and decoUimation observed 
in other sources (e.g. 3C31). These two sources arc repre- 
sentative examples of two kinds, which are classified into 
Fanaroff-Riley type I (FRI) and Fanaroff-Riley type II 
(FRIT) radio sources (FanarofF & Riley 1974). The mor- 
phological dichotomy of FRI and FRII sources may be 
related to the stability properties of relativistic jets with 
different kinetic powers. 

This complex situation motivated us to study the in- 
teraction of relativistic jets with their ambient medium 
and more specifically the KH instability in detail, by ap- 
plying a linear stability analysis along with numerical sim- 
ulations. 

The linear analysis of KH instability in relativistic jets 
started with the work of Turland & Scheuer (1976) and 
Blandford & Pringle (1976) who derived and solved a dis- 
persion relation for a single plane boundary between the 
relativistic flow and the ambient medium. Next, Ferrari, 
Trussoni & Zaninetti (1978) and Hardee (1979) examined 
properties of KH instability in relativistic cylindrical jets 
by following the derivation of the dispersion relation in 
the nonrelativistic case, in the vortex sheet approxima- 
tion, done by Gill (1965). They numerically solved the 
dispersion relation, found unstable KH modes and clas- 
sified them into the fundamental (surface) and reflection 
(body) family of modes. The classification is related to the 
number of nodes, across the jet, of sound waves reflecting 
in between jet boundaries. An internal wave pattern is 
formed by the composition of oblique waves, for which 
the jet interior is a resonant cavity. The physical mean- 
ing of KH instability in supersonic jets has been discussed 
by Payne & Cohn (1985), who have shown that the pres- 
ence of instability is associated with the overreflection of 
soundwaves (the modulus of reflection coefficient is larger 
than 1) on the sheared jet boundaries. 

Subsequent studies include effects of magnetic field 
(Ferrari, Trussoni & Zaninetti 1981), the effects of the 
shear layer, replacing the vortex sheet in the nonrelativis- 
tic planar case (Ferrari, Massaglia & Trussoni 1981; Ray 
1982: Choudhury & Lovelace 1984), and conical and cylin- 
drical jet geometry (Birkinshaw 1984; Hardee 1984, 1986, 
1987). The effects of cylindrical a shear layer have been 
examined by Birkinshaw (1991) in the nonrelativistic case, 
attempted by Urpin (2002) in the relativistic case and the 
presence of multiple components (jet + sheath + ambient 
medium) was investigated by Hanasz & Sol (1996, 1998). 
Other authors have investigated current-driven modes in 
magnetized jets (Appl & Camenzind 1992; Appl 1996) in 
addition to KH modes. 

An extension of the linear stability analysis in the rel- 
ativistic case to the weakly nonlinear regime has been per- 
formed by Hanasz (1995) and led to the conclusion that 
Kelvin- Helmholtz instability saturates at finite amplitudes 
due to various nonlinear effects. An explanation of the na- 



ture of the mentioned nonlinearities has been proposed by 
Hanasz (1997). The most significant effect results from the 

relativistic character of the jet flow, namely from the fact 
that the velocity perturbation cannot exceed the speed of 
light. 

A more recent approach is to perform a linear stabil- 
ity analysis in parallel with numerical simulations and to 
compare the results of both methods in the linear regime 
and then to follow the nonlinear evolution of the KH in- 
stability resulting from numerical simulations. Hardee & 
Norman (1988) and Norman & Hardee (1988) have made 
such a study for nonrelativistic jets in the spatial approach 
and Bodo et al. (1994) in the temporal approach. In the 
relativistic case this kind of approach was applied for the 
first time by Hardee et al. (1998) in the case of axisym- 
metric, cylindrical jets and then extended to the 3D case 
by Hardee et al. (2001). 

Similarly to the linear stability analysis, the numeri- 
cal simulations of jet evolution can be performed following 
both the spatial and the temporal approach, depending on 
the particular choice of initial and boundary conditions. 
In the temporal approach one considers a short slice of 
jet limited by periodic boundaries along the jet axis, and 
adds some specific perturbation, eg. an eigenmode result- 
ing from the linear stability analysis. Due to the periodic 
boundary conditions the growing perturbations can only 
be composed of modes having a wavelength equal to the 
length of the computational box and/or its integer frac- 
tions (Bodo et al. 1994, 1995, 1998). Whereas the spatial 
approach appears more appropriate to analyze the global 
dynamics and morphology of the whole jet, the temporal 
approach is suitable for the comparison between the nu- 
merical results and analytical studies of the jet stability 
because, due to the fact that only part of the jet is sim- 
ulated, a high effective numerical resolution is achievable 
with limited computer resources. 

Numerical simulations (Marti et al. 1997; Hardee et 
al. 1998, Rosen et al. 1999) demonstrate that jets with 
high Lorcntz factors and high internal energy htq influ- 
enced very weakly by the Kelvin-Hclmlioltz instability. 
Moreover, Hardee et al. (1998), Rosen et al. (1999) note 
that contrary to the cases with lower Lorcntz factors and 
lower internal energies, the relativistically hot and high 
Lorentz factor jets do not develop modes of KH instabil- 
ity predicted by the linear theory. They interpret this fact 
as the result of a lack of appropriate perturbations gen- 
erating the instability in the system. In the limit of high 
internal energies of the jet matter the Kelvin-Helmholtz 
instability is expected to develop with the highest growth 
rate. 

In this paper we focus on the role of internal energy and 
the Lorentz factor for the stability of relativistic planar 
two-dimensional jets in the temporal approach. We chose 
this simple configuration in order to enhance our ability 
to understand the complex nonlinear evolution of the un- 
stable KH modes in relativistic jets. For this purpose, we 
perform a linear analysis of the KH instability and con- 
struct perturbations analytically: the eigenmodes of the 
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linear problem representing the most unstable first reflec- 
tion mode for each set of jet parameters. These modes 
are then incorporated as initial states of numerical simu- 
lations. The choice of planar two-dimensional jets allows 
us to study symmetric and antisymmetric modes of KH 
instabihty which resemble, respectively, pinch and heli- 
cal modes of cylindrical jets. Only the flutting modes of 
cylindrical jets, corresponding to azimuthal wavenumbers 
larger than 1 , have no counterparts among the eigenmodes 
of planar two-dimensional jets. In the present paper we in- 
vestigate only the symmetric modes. 

Our aim is to investigate physical details of the process 
of transition of KH instability modes from the linear to the 
nonlinear regime. In the next step (Paper II) we shall in- 
vestigate the relations of long-term nonlinear evolution of 
KH instability (including the mixing phase and the prop- 
erties of the evolved flow at the end of our simulations) on 
the initial jet parameters. 

The paper is organized as follows. In Section |21 we de- 
scribe our method, i.e. perform the linear analysis of the 
KH instability for the relativistic planar case, then define 
perturbations and describe numerical algorithm and other 
details of numerical simulations. In Section |3| we present 
and discuss results of the both linear stability analysis and 
the numerical simulations. In Scction^jwc summarize and 
conclude our paper. 

2. Method 

Following the standard procedure (see eg. Gill 1965; 
Ferrari, Trussoni & Zaninetti 1978; Hardee 1979) we derive 
the dispersion relation for the Kelvin-Helmholtz modes. 
We focus on the simplest geometrical configuration of two- 
dimensional planar relativistic jets and apply the temporal 
stability analysis. 



2.1. Linear stability analysis 

The full set of equations describing the current problem 
consists of the set of relativistic equations of hydrodynam- 
ics for a perfect fluid, (e.g. Ferrari, Trussoni & Zaninetti 
1978) 
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dv 



_ V dp 



= 0, 



fdp 



'dt 



+ V- (71;) 



and an (adiabatic) equation of state 



PPq = const. 



(1) 



(2) 



(3) 



In the preceding equations, c is the speed of light, po is the 
particle rest mass density (i.e., pq = mn, where m is the 
particle rest mass and n the number density in the fluid 
rest frame), p stands for the relativistic density which is 
related to the particle rest mass density and the specific 
internal energy, e, by p = po{l -\- e/c^). The enthalpy is 



defined as w — p + pjc^, the sound speed is given by 
Cj = Tp/w, where F is the adiabatic index. The relation 
between pressure and the specific internal energy is p — 
(r — l)epo- The velocity of the fluid is represented by v 
and 7 is the corresponding Lorentz factor. 

The assumed geometry of the jet considered in the 
forthcoming linear stability analysis and the numerical 
simulations is sketched in Fig. ^ First of all the consid- 
ered jet is 2D-planar and symmetric with respect to the 
X = plane. The flow in the jet moves in the positive 
z direction and its matter forms a contact discontinuity 
(a vortex sheet) with the matter of external medium at 
X = —Rj and x = Ry From now on all quantities repre- 
senting the jet will be assigned with the 'j' subscript and 
the quantities representing the ambient medium will be 
assigned with 'a'. 

The following matching conditions are imposed on the 
interface between the jet and the ambient medium 



Pa 

h„ ■ 



Pj for 
hj for 



Rj 
Rj. 



(4) 
(5) 



The matching conditions express the assumption of equal- 
ity of the jet and ambient pressures and the equality of 
transversal displacements of jet {hj) and ambient (ha) 
fluid elements adjacent to the jet boundary (at — Rj). 

In addition, the Sommerfeld radiation condition (ex- 
pressing the disappearance of perturbations at the infin- 
ity) will be applied for linear perturbations. 



2.1.1. Equilibrium state for the linear stability analysis 

We assume that the initial state is an equilibrium configu- 
ration. The initial equilibrium can be described by the fol- 
lowing set of independent parameters: the Lorentz factor 
corresponding to the unperturbed longitudinal jet veloc- 
ity, Vj, 7 = (1 — u|/c^)~^/^, the particle rest mass density 
of the jet poj (the particle rest mass density of the ambient 
medium is normalized to unity: poa = 1) and the specific 
internal energy of the jet Sj. The ambient medium is as- 
sumed to be at rest {va — 0). 

The other dependent parameters describing the 
equilibrium state are: the internal jet Mach number 
Mj = ''^j/csj corresponding to the initial jet lon- 
gitudinal velocity, the relativistic rest mass density 
Poj{l+ej/c^) 



contrast v ■ 



POa(l + £a/c^) 



(or, equivalently, the en- 



thalpy contrast rj — 



po,(l + (r,-l)e,/cV(l + £,/c2)) 



POa (1 + (Ta - l)ea/cV(l + Ea/C^)) ' 

and the specific internal energy of the ambient medium 
Ea — 777^ r-^Ej, which is related to the specific inter- 

[Ta — l)P0a 

nal energy of the jet through the pressure balance condi- 
tion. 
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(extended grid to 100-200 Rj ) 



3 Rj uniform grid (400 cells/Rj) 



(1) 



z (flow direction) 



XRj (16cells/R|) 



2Ri 



Rj 



jet/ambient medium 
boundary layer 

(1) 



-jet axis (reflection boundary conditions) 



(1) periodic boundary conditions 



Fig. 1. The geometry of the flow considered in the linear stability analysis and the numerical simulations, including a 
description of the boundary conditions. 



2.1.2. Dispersion relation 

The first step towards the dispersion relation is to reduce 
the equations iQ) - |(SJ| to the dimensionless form through 
the normalization of spatial coordinates to the jet radius 
Rj, velocities to the sound speed of the jet material Cgj, 
time to the dynamical time Rj/cgj and pressure to the 
equilibrium pressure. 

The next step is to decompose each dependent quan- 
tity into the equilibrium value and the linear perturba- 
tion. After the reduction of equations to the dimensionless 
form and substitution of the perturbed quantities in equa- 
tions (0-©, we obtain the following set of dimensionless 
linearized equations. In the following the dimensionless 
quantities will be assigned the same symbols as the pre- 
vious dimensional ones. The subscript '1' stands for the 
linear perturbation of the corresponding variable. For the 
jet medium we get 

Vj dp J I 



r.-7'f^ + K-V)t;,i)+Vp,i 



V dt 



dt 



dp 



dt 



r,V-t;,i 



= 0, (6) 



(7) 



^,2r "'J 



dt 



where Vj is the initial (unperturbed) jet velocity in units 
of jet internal sound speed, c is the speed of light in units 



of the sound speed and the normalized pressure Pq = 1 \s 
omitted in the equations. For the ambient medium we get 



Fj dvgi 
T] dt 

dpal 

dt 



(8) 



(9) 



The linearized matching conditions 10} and Q at the jet 
interface read 



Pal = Pji for \x\ = 1 
hai = hji for |a;| = 1, 



(10) 
(11) 



where the displacements of fluid elements adjacent to 
the contact interface in the linear regime are related to 
transversal velocities by the following formulae 



V-jxl 



d 



d 



dt ' dz 



Vaxl — Trrhai- 
dt 



(12) 
(13) 



The following wave equations can be derived from the 
equations © - respectively for the jet 
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Vaxi = -iu;h.,i. 



(22) 



d_ 
di 



Pji 



(14) 



-7 



d ^Vj d —0 
dz dt ' ^"'^ 



and for the ambient medium 

d'^Pal rjTa ( d'^Pal , d^Pal 



a<2 



r, 



da 



dz 



0. 



(15) 



It is apparent that Eqn. H14|l and (|15() describe propaga- 
tion of obhque sound waves in the jet and ambient medium 
respectively. 

The next step is to substitute perturbations of the form 
Sji = [S^F+{x) + SjFr{x)] expi^z - cut) + c.c. (16) 

with Fp{x) — exp(±i /cjxa;) to describe waves propagat- 
ing in positive and negative x-directions in the beam. 
Here we use fc|| and ka^ji_ for longitudinal and transverse 
wavenumbers. Perturbations in the external medium are 
of the form 



5al = SaFai^) "^^P i{k\\Z — Ljt) + C.C. 



(17) 



where Fi^{x) = exp(ika±x) for |a;| > 1. Assuming that 
Re(fcQ^) > 0, only outgoing waves are present in the am- 
bient medium. After the substitution of the explicit forms 
of pressure perturbations to equations H14|) and IjlSfl we ob- 
tain the following expressions for the perpendicular com- 
ponents of wavevectors 



kaA 



- J , ,2 



1/2 



kj± = (u;' -k'/) 



/2.1/2 



(18) 
(19) 



which are standard relations for linear sound waves in 
both media (note that rjTa/Tj is the squared sound speed 
of ambient medium in units of the sound speed of jet), 
to' = — Vjk\^) and fcj| = 7(fc|| — ^uj) are frequency 
and wavenumber of the internal sound wave in the ref- 
erence frame comoving with the jet. The wavevectors 
kj — (fcjj^, /cjIi) and ka = (fca±,fcajt) determine the direc- 
tion of propagation of sound waves in the jet and ambient 
medium respectively. Vanishing ofkjj_, for instance, would 
mean that the jet internal sound waves move in the axial 
direction. In cases kj± ^ and/or ka± ^ the propaga- 
tion of the sound waves is oblique with respect to the jet 
axis. 

The corresponding perturbation of the jet surface can 
be written as 



hai = hji = Ziexp i(fc||Z — cut) + c.c. 



(20) 



which after the substitution to equations H12(l and (|13|l 
reads 



jxi = -i{uj - Vjk\i) hji, 



(21) 



With the aid of the above perturbations the whole sys- 
tem of partial differential equations is reduced to a set of 
homogeneous linear algebraic equations. The dispersion 
relation appears as solvability condition of the mentioned 
set of equations, namely it arises from equating the deter- 
minant of the linear problem to zero. Within the present 
setup, the dispersion relation for the Kelvin-Helmholtz in- 
stability in supersonic, relativistic, two-dimensional slab 
jets can be written as 



1 (C.'^ - fc'jj)l/2 



1/2 



(23) 



fcoth i{uj'^ - A:'jj)i/2 for s = 1 
,th^(cJ'^-fc'jj)l/2 for s = -l 

where s = ±1 is the symmetry of perturbation. In the 
present paper we focus on symmetric perturbations (s = 
l). 

The above dispersion relation is solved numerically 
with the aid of the Newton-Raphson method (see, e.g., 
Press et al. 1992). In this way we obtain the complex 
frequency as a function of the parallel component of the 
wavenumber /c|| , since we choose the temporal stability 
analysis for the present investigations. 



2.1.3. Eigenstates of the system 

Once the solutions of the dispersion relation Lu{k) are 
found, the derivation of eigenstates of the system is 
straightforward. By utilizing the set of relativistic equa- 
tions ® - © one can relate perturbations of gas density 
and velocity to the perturbation of pressure. 

First, we choose the amplitude of the pressure wave in 
the jet, p^, as a free constant and then relate the other 
constants to its value. The corresponding '— ' amplitudes 
are computed by multiplying these values by s in case of 
pressure and longitudinal velocity and by — s in the case 
of transversal velocity. 

The amplitude of pressure wave in the ambient 
medium, is found from the pressure matching con- 
dition at a; = 1: 



p+ exp(i kai_) = pf exp(i kj±) + sp. exp(-i kj 



(24) 



The remaining amplitudes of velocity perturbations are 



+ _ r]_kai_ + 

Pa ' 



IP J ' 



(25) 
(26) 
(27) 
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where 77 is the ratio of enthalpies as defined earlier. 

Regarding specific internal energy and rest-mass den- 
sity, perturbations are calculated as follows: 



-Pajl 



Pa,jl 



1 



Pa,j + l/c^ Ta, 



(29) 



(30) 



where paj refer to the relativistic density in ambient and 
jet media in dimensionless units, respectively, and Paji is 
the corresponding perturbation. 

2.1.4. Lorentz transformation of velocity perturbations 

In the forthcoming discussion of results of numerical simu- 
lations we shall analyze the velocity components in the ref- 
erence frame comoving with the jet. The standard Lorentz 
transformation formulae for the velocity components are 



1 - VzVj/c^ 



^a: -1 / 2 

1 - VzVj/c^ 



(31) 



Assuming that the Lorentz transformation is applied for 
linear perturbations of sufficiently small amplitude, these 
transformation rules can be approximated by the following 
formulae 



/ 2 



41 



(32) 



where Vxl^ Vzi are the velocity perturbations in the ref- 
erence frame of the ambient medium and 



^xlJ "zl 



are 



the velocity perturbations in the jet reference frame. 
These transformation rules tell us that the longitudinal 
and transversal velocity perturbations in the jet refer- 
ence frame are larger by factors of 7^ and 7 respectively 
than the corresponding components in the rest (ambient 
medium) frame. 

2.2. Numerical simulations 

Now we proceed by performing numerical simulations of 
jet models to study the growth of unstable modes through 
the linear and non-linear phases. To this aim, we adopt 
as an equilibrium initial state the one described in the 
former section. In order to avoid the growth of random 
perturbations with wavelengths of the order of the cell 
size, we replace the transverse discontinuous profiles of 
equilibrium quantities by smooth profiles of the form 



Vz{x) 



cosh(a;™) 



po{x) = poa 



POa - POj 

cosh(a;™) 



(33) 



(34) 



(see Bodo et al. 1994), where m is the steepness parameter 
for the shear layer. Typically we used m — AO, for which 
the shear layer has a width ~ O.lRj (see Appendix A). 

On the other hand the linear solutions derived in 
Section l2. 1 .31 correspond to an idealized case with a con- 
tact discontinuity at the jet interface. Finding of analogous 
solutions for the corresponding sheared boundary problem 
is more difficult from the numerical point of view (see Ray 
1982, and Choudhury & Lovelace 1984), therefore we de- 
cided to implement corresponding smoothed eigenstates 
obtained for the vortex sheet limit as approximate solu- 
tions for the case of sheared boundary. The validity of 
this procedure is proven a posteriori by the convergence 
of theoretical growth rates and the ones determined for 
simulated KH modes. Since, in case of the vortex sheet 
approximation, all first order dynamical variables (except 
pressure and the transversal displacement of fiuid element 
adjacent to the jet boundary) are discontinuous, we apply 
the following smoothing of the equilibrium quantities and 
linear perturbations at the thin layer surrounding the jet 
boundary 



5l = Sal - 



Sal - Sji 

cosh(a;™) ' 



(35) 



where S stands for any perturbed variable. Such smoothed 
eigenstates are subsequently implemented as small- 
amplitude perturbations in the initial equilibrium model 
of the simulations. 

Our first aim is to reproduce the linear evolution of 
unstable modes by means of hydrodynamical simulations. 
For this reason and in order to make the picture of evolu- 
tion of the KH instability as clear and simple as possible, 
we perturb the equilibrium with the most unstable first 
reflection mode, which is smoothed at the jet interface ac- 
cording to formula (|35|l . The length of the computational 
grid is in each simulation equal to the wavelength of the 
mentioned mode. Defining the initial state in such a way, 
we know precisely which growth rate to expect in numer- 
ical simulations. The difference between the expected and 
computed growth rates is a measure of the performance 
of the code in describing the linear regime and gives us 
confidence on the numerical results from the non-linear 
evolution. 

The simulations have been performed using a high- 
resolution shock capturing code to solve the equations 
of relativistic hydrodynamics in Cartesian planar coor- 
dinates. The code is an upgrade of that used to study 
large scale (Marti et al. 1997) as well as compact relativis- 
tic jets (Gomez et al. 1997). The equations of relativistic 
hydrodynamics are written in conservation form and the 
time variation of the proper rest mass, momentum and 
total energy within numerical cells are calculated through 
the fluxes across the corresponding cell interfaces. Fluxes 
are calculated with an approximate Riemann solver that 
uses the complete characteristic information contained in 
the Riemann problems between adjacent cells (Donat & 
Marquina 1996). It is based on the spectral decomposition 
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of the Jacobian matrices of the relativistic system of equa- 
tions derived in Font et al. (1994) and uses analytical ex- 
pressions for the left eigenvectors (Donat et al. 1998). The 
spatial accuracy of the algorithm is improved up to third 
order by means of a conservative monotonic parabolic re- 
construction of the pressure, proper rest-mass density and 
the spatial components of the fluid four-velocity, as in 
GENESIS (see Aloy et al. 1999, and references therein). 
Integration in time is done simultaneously in both spa- 
tial directions using a total-variation-diminishing (TVD) 
Runge-Kutta scheme of high order (Shu & Osher 1988). 
See Marti et al. (1997) for the differential equations as 
well as their finite-difference form, and a description of 
exhaustive testing of the hydrodynamical code. Besides 
relativistic density, momentum and energy, the code also 
evolves a passive scalar representing the jet mass fraction. 
This allows us to distinguish between ambient and jet mat- 
ter helping us to characterize processes like jet/ambient 
mixing or momentum exchange. 

The parameters of the simulations presented in this 
paper are listed in Tabled Values of the parameters were 
chosen to be close to those used in some simulations by 
Hardee et al. (1998) and Rosen et al. (1999) and are chosen 
to span a wide range in thermodynamical properties as 
well as beam flow Lorentz factors. In all simulations, the 
density in the jet and ambient gases are poj =0.1, poa = 1 
respectively and the adiabatic exponent Tj_a — 4/3. 

Since the internal rest mass density is flxed, there are 
two free parameters characterizing the jet equilibrium: the 
Lorentz factor and the jet speciflc internal energy dis- 
played in columns 2 and 3 of Tabled Models whose names 
start with the same letter have the same thermodynamical 
properties. Beam (and ambient) speciflc internal energies 
grow from models A to D. Three different values of the 
beam flow Lorentz factor have been considered for models 
B, C and D. The other dependent parameters mentioned 
in Section r2. 1.11 are displayed in columns 4-10 of Tabled 
Note that given our choice of poj , the ambient media as- 
sociated to hotter models are also hotter. The next three 
columns show the longitudinal wavenumber together with 
oscillation frequency and the growth rate of the most un- 
stable reflection mode. The following three columns dis- 
play the same quantities in the jet reference frame. Next 
two columns show the perpendicular wavenumbers of lin- 
ear sound waves in jet and ambient medium respectively. 
The last column shows the linear growth rate of KH in- 
stability in the jet reference frame expressed in dynamical 
time units, i.e. in which time is scaled to Rj/Cgj- All other 
quantities in the table, are expressed in units of the ambi- 
ent density, poa, the speed of light, c, and the jet radius, 
R,. 

As it is apparent in Table ^ growth rates tend to in- 
crease with the specific internal energy of the beam and to 
decrease with the beam flow Lorentz factor. Note that, in 
the jet reference frame, models with the same thermody- 
namical properties tend to have (within « 10 %) the same 
growth rates. Note also that in the jet reference frame and 



in dynamical time units, all the models have comparable 
(within a factor 1.5) linear growth rates. 

The initial numerical setup consists of a steady two- 
dimensional slab jet model (see Fig. As stated above, 
a thin shear layer between the ambient medium and the 
jet is used instead of the vortex sheet. Due to symmetry 
properties, only half of the jet {x > 0) has to be com- 
puted. Reflecting boundary conditions are imposed on the 
symmetry plane of the flow, whereas periodical conditions 
are settled on both upstream and downstream boundaries. 
The numerical grid covers a physical domain of one wave- 
length along the jet (3 to 20 Rf, see Table [TJ and 100 Rj 
across (200 Rj in the case of models D). The size of the 
transversal grid is chosen to prevent loses of mass, momen- 
tum and energy through the boundaries. 400 numerical 
zones per beam radii are used in the transverse direction 
across the flrst 3 Rj. From this point up to the end of 
the grid, 100 (200, in case of models D) extra numerical 
zones growing geometrically have been added. The width 
growth factor between contiguous zones is approximately 
1.08 for models A, B and C and 1.04, for models D. Along 
the jet, a resolution of 16 zones per beam radius has been 
used. The agreement of the linear stability analysis with 
the results of numerical simulations depends on the grid 
resolution (and the width of the initial shear layer). The 
applied resolution of 400 x 16 grid zones per beam radius 
is chosen on the base of several tests, which are presented 
in detail in Appendix A. 

The steady model is then perturbed according to the 
selected mode, with an absolute value of the pressure per- 
turbation amplitude inside the beam of = 10~^. This 
means that those models with the smallest pressure, like 
model A, have relative perturbations in pressure three or- 
ders of magnitude larger than those with the highest pres- 
sures, D. However this difference seems not to affect the 
linear and postlinear evolution. 

3. Results 

Following the behaviour of simulated models we flnd that 
the evolution of the perturbations can be divided into a 
linear phase, a saturation phase and mixing phase. This 
section is devoted to describing the properties of the linear 
and saturation phases in our models focusing on the influ- 
ence of the basic parameters. The remaining fully nonlin- 
ear evolution will be discussed in Paper II. Our description 
shares many points with that of Bodo et al. (1994) for the 
case of classical jets. 

In order to illustrate the growth of perturbations and 
determine the duration of the linear and saturation phases 
in our simulations, we plot in Fig.Elthe amplitudes of the 
perturbations of the longitudinal and transversal veloci- 
ties, inside the jet and in the jet reference frame, together 
with the pressure oscillation amplitude. We plot also the 
growth of the imposed eigenmodes resulting from the lin- 
ear stability analysis. Both the velocity perturbations are 
transformed from the ambient rest frame to the unper- 
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THyn 



Model 


7 


J 


A05 


5 


0.08 


B05 


5 


0.42 


C05 


5 


6.14 


DOS 


5 


60.0 


BIO 


10 


0.42 


CIO 


10 


6.14 


DIO 


10 


60.0 


B20 


20 


0.42 


C20 


20 


6.14 


D20 


20 


60.0 



p 



Mi 



fell 



0.008 
0.042 
0.614 
6.000 
0.042 
0.614 
6.000 
0.042 
0.614 
6.000 



0.18 
0.35 
0.55 
0.57 
0.35 
0.55 
0.57 
0.35 
0.55 
0.57 



0.059 
0.133 
0.387 
0.544 
0.133 
0.387 
0.544 
0.133 
0.387 
0.544 



0.0027 
0.014 
0.205 
2.000 
0.014 
0.205 
2.000 
0.014 
0.205 
2.000 



0.11 0.11 
0.14 0.15 
0.44 0.51 
0.87 0.90 
0.14 0.15 
0.44 0.51 
0.87 0.90 
0.14 0.15 
0.44 0.51 
0.87 0.90 



5.47 
2.83 
1.80 
1.71 
2.88 
1.83 
1.73 
2.89 
1.83 
1.74 



0.30 
0.69 
2.00 
2.63 
0.50 
1.91 
2.00 
0.46 
1.44 
2.00 



0.20 0.026 

0.49 0.055 

1.60 0.114 

2.18 0.132 

0.41 0.031 

1.72 0.055 

1.81 0.063 

0.39 0.014 

1.37 0.027 

1.91 0.029 



1.32 
2.62 
5.73 
7.02 
3.59 
9.77 
9.67 
6.51 
13.89 
18.11 



7.20 0.13 

7.32 0.28 

9.98 0.57 

11.56 0.66 

10.28 0.31 

17.67 0.55 

16.58 0.63 

18.76 0.28 

25.38 0.54 

31.43 0.58 



7.08 
6.84 
8.17 
9.18 
9.64 
14.72 
13.47 
17.60 
21.24 
25.68 



0.53 
1.08 
1.07 
0.24 
0.94 
1.49 
0.20 
0.90 
1.28 
0.31 



0.73 
0.79 
1.05 
1.15 
0.90 
1.01 
1.10 
0.81 
0.99 
1.01 



Table 1. Equilibrium parameters of different simulated jet models along with solutions of the dispersion relation 
H24(l . corresponding to fastest growing first reflection mode, taken as input parameters for numerical simulations. The 
meaning of the symbols is as follows: 7: jet flow Lorentz factor; e: specific internal energy; Cg'- sound speed; pp: pressure; 
v. jet-to-ambient relativistic rest mass density contrast; rj: jet-to-ambient enthalpy contrast; Mj: internal jet Mach 
number; k\\ _i: longitudinal/transverse wavenumbers; uJr,i- real/imaginary part of the wave frequency. Labels a and 
j refer to ambient medium and jet, respectively. The primes are used to assign wavenumber and complex frequency 
in the reference frame comoving with jet. The last column shows the linear growth rate of KH instability in the jet 
reference frame expressed in dynamical time units, i.e. in which time is scaled to Rj /csj. All the quantities in the table, 
except the last column, are expressed in units of the ambient density, poa, the speed of light, c, and the jet radius, Rj. 



turbed jet rest frame using the Lorentz transformation 
rules for velocity components, given by Eqn. H31|) . 

We define the characteristic times in the following way. 
During the linear phase the ratios of the exponentially 
growing amplitudes of pressure, longitudinal velocity and 
transversal velocity remain constant by definition. We de- 
fine the end of the linear phase (tun) as the moment at 
which one of quantities starts to depart from the initial 
exponential growth. Within the set of our models the first 
quantity to break the linear behaviour is the longitudinal 
velocity perturbation. 

Later on the transversal velocity saturates, i.e. stops 
growing at the saturation time (tsat)- We call the period 
between and isat the saturation phase. We find also 
that the pressure perturbation amplitude reaches a max- 
imum. This moment is denoted by ipoak- In Paper II we 
will see that this peak announces the entering of the fully 
nonlinear regime. The choice of tun, ^sat, ipcak has been 
illustrated in Fig. [3 (top panel). 

Table |21 collects times of the linear and saturation 
phases in the different models (columns 2-4) . The last col- 
umn shows the saturation time in dynamical units and 
in the jet reference frame. The change of reference frame 
eliminates the effects coming from the jet flow Lorentz fac- 
tor that stretches out the rhythm of evolution in the am- 
bient rest frame. Dynamical time units are adapted to the 
characteristic time of evolution of each model. Hence this 
change of units and reference frame allows us to compare 
the relevant scales of evolution of all the models directly. 

3.1. Linear phase 

In the linear phase the ratios of oscillation amplitudes of 
different dynamical variables (density, pressure and veloc- 
ity components) are constant. We have used this property 



Model 


ilin 


^sat 


^pcak 


,/dyn 

^ sat 


A05 


180 


380 


380 


13.69 


B05 


125 


200 


205 


13.99 


C05 


100 


125 


130 


13.74 


D05 


105 


120 


130 


13.71 


BIO 


235 


380 


385 


13.29 


CIO 


210 


245 


250 


13.46 


DIO 


180 


225 


225 


12.86 


B20 


450 


760 


780 


13.29 


C20 


270 


645 


775 


17.72 


D20 


350 


480 


500 


13.71 



Table 2. Times for the different phases in the evolution 
of the perturbed jet models. In columns 2-4 we show: <iin - 
end of linear phase (the ratios of amplitudes of the differ- 
ent quantities are not constant any longer), ^sat - end of 
saturation phase (the amplitude of the transverse speed 
perturbation reaches its maximum), tpcak - the peak in 
the amplitude of the pressure perturbation is reached, (see 
Fig. 121 . Note that, as a general trend, < ts&t < ^poak 
. In columns 5 we show the saturation time in dynamical 
time units and in the jet reference frame. 



to identify the end of the linear phase with the time tun at 
which one of the variables deviates from linear growth in a 
systematic way. Note that our definition is different from 
that of Bodo et al. (1994) who associate the end of the 
linear phase with the formation of the first shocks inside 
the jet. 

We note that in our simulations, during the linear 
phase, the growth of perturbations of each dynamical vari- 
able follows the predictions of the linear stability analysis 
with relatively good accuracy. On average the growth rates 
measured for all our numerical experiments are about 20 
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% smaller than corresponding growth rate resulting from 
the linear stability analysis. This small discrepancy may 
be partially a result of the application of the shearing layer 
in simulations and the vortex sheet approximation in the 
linear stability analysis, and also to a lack of transversal 
resolution in the numerical simulations (see Appendix A). 

In all cases the amplitude of the longitudinal velocity 
oscillation is the first quantity to stop growing. The rea- 
son for the limitation of the velocity oscillations is obvious: 
the oscillations of velocity components (corresponding to 
sound waves propagating in the jet interior) cannot ex- 
ceed the speed of light. This limitation (specific to rela- 
tivistic dynamics) is easily noticeable in the jet reference 
frame, but it is obscured by the Lorentz factor (in the first 
and second power) in the reference frame of the ambient 
medium, as it follows from the formulae H32|) . 

As seen in Fig. [21 model C20 evolves distinctly from the 
other models: at t = 270, the linear phase ends before the 
saturation of longitudinal speed. See the middle panels 
of Fig. [71 (especially the pressure panel) which could be 
responsible for the peculiar evolution of model C20. Note 
also that model C20 is the only one that has a theoretical 
perturbation growth rate smaller than the numerical one 
(probably associated with the excited short wavelength 
mode). 

At the end of the linear phase, the values of quanti- 
ties like density, pressure and flow Lorentz factor are still 
very close to the corresponding background values (the 
perturbation in pressure is between 10 % and 50% of the 
background pressure in all the models). Figure [3 shows 
the distribution of pressure at the end of the linear phase 
for two representative models BOS and D05. 

3.2. Saturation phase 

As it has been mentioned, the end of the linear phase co- 
incides with the limitation of the longitudinal oscillation 
velocity. At time tu^ the transversal velocity component 
is smaller than the speed of light, by an order of magni- 
tude approximately, so the transversal velocity perturba- 
tion has still room to grow. We defined tsat as the time 
corresponding to the saturation of the transversal velocity 
and saturation phase as the period between tn^ and tgat- 
We find that the longitudinal velocity perturbation 
amplitude reaches almost the speed of light (more than 
0.9c) while the transversal perturbation amplitude stops 
its growth at the level of approximately 0.5c for all pre- 
sented simulations. This can be explained in the following 
way. We remind that the eigenmodes are built of oblique 
sound waves overreflected at the jet/ambient medium in- 
terface. This means that the amplitude of the reflected 
wave is larger than the amplitude of the incident one. 
Sound waves are longitudinal waves, so the gas oscillation 
velocity can be split into the longitudinal and transver- 
sal components separately for the incident and reflected 
waves. Let us consider locally the incident wave and the re- 
flected wave in the jet medium, in a flxed point close to the 



jet boundary. Then the longitudinal velocity components 
of the incident and reflected waves sum in phase, while 
the transversal ones sum in counterphase. This means that 
while the amplitude of the total longitudinal velocity oscil- 
lation component grows and approaches the speed of light, 
the oscillation amplitude of the total transversal compo- 
nent is a difference of the two velocities, each smaller than 
the speed of light. Therefore the oscillation amplitude of 
the total transversal component has to take values signif- 
icantly smaller than the longitudinal one. 

The duration of the saturation phase depends on the 
Lorentz factor and the speciflc internal energy with a ten- 
dency to increase with the former and to decrease with the 
latter (exception made of models CIO or DIO and C20). 
It ranges between a few tens of (absolute) time units to a 
few hundreds. 

We note that the saturation times expressed in dy- 
namical time units and in the jet reference frame (see 
Table [JJ are almost equal for all models. This similar- 
ity can be explained by the following argument. First, 
the amplitude of pressure perturbations pf = 10~^ is the 
same for all simulations. Second, the amplitudes of pres- 
sure and transversal velocity perturbations are related by 
formula (|26|l . The perturbation of perpendicular velocity 
is then transformed to the jet reference frame following 
(I32|l . Following the results of numerical experiments, the 
transversal velocity grows with the linear growth rate un- 
til the transversal velocity perturbation reaches the upper 
limit, i.e. exp(a;itsat) — 0.5c. If we express tgat in dy- 
namical time units then in the jet reference frame 



7 



CsjUJi 



-if 



(36) 



Since the term uj' /kj± under the logarithm is close to the 
jet sound speed and is varying only slightly in between 
models, while the other factors under the logarithm are 
constants, the saturation time is proportional to the in- 
verse growth rate. From the linear stability analysis, the 
latter one measured in the jet reference frame and ex- 
pressed in dynamical time units is almost equal for all 
the models (see last column in Table QJ. Substitution of 
numbers into the formula (|36() provides t'^^^ in the range 
10-^ 15, which is consistent with the value 13.5 in Table[21 
However, the remarkable convergence of saturation times 
for all the models is difficult to explain in view of the 
additional randomness in the evolution of perturbations, 
which is apparent in Fig. [2 

During the saturation phase, the jet inflates (Bodo et 
al. 1994 called this phase the expansion phase) and de- 
forms due to transversal oscillations. On the other hand, 
the saturation time igat coincides (within a few time units) 
with the appearance of an absolute maximum in the pres- 
sure distribution (at tpoak) at the jet boundary, and with 
the start of the mixing phase (to be studied in Paper II). 
Figure EE display snapshots of several quantities close to 
the end of the saturation phase before the distortion of 
the jet boundary. 
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Fig. 2. Evolution of the relative amplitudes of perturbations. Dotted line: pressure perturbation ((p„ 



Po)/po)- 



Dashed line: longitudinal velocity perturbation in the jet reference frame (0.5 — v'^ mm))- Dash-dotted line: 

perpendicular velocity perturbation in the jet reference frame (0.5 {v'^ ^^j. — v'^ min))- The search for maximum {{pmax, 
Wx,max-i Wz,max) ^'^'^ minimum (w^^mm: Wz,min) values have been restricted to those numerical zones with jet mass 
fraction larger than 0.5. Solid line: linear analysis prediction for the growth of perturbation. Note that abscissae in 
the last column plots extend up to t = 1000i?j/c and that ordinate values adapt to fit the scale of each plot. Arrows 
in the plot of Model A05 point to specific stages of evolution used to define tun, isat and tpeak (see text). 



The structure of perturbations at the end of the sat- 
uration phase is quite similar in all models. The longitu- 
dinal wavelength of perturbations is given by the linear 
stability analysis and it is constant because of the fixed 
length of the computational domain. Along with the lon- 
gitudinal wavelength, the opening angle of oblique waves, 
far from the jet symmetry plane, is given by the linear 
analysis. The opening angle of oblique waves is enlarged 
in the vicinity of the jet interface. Closer to the jet sym- 
metry plane the perturbations are already in the nonlinear 
phase, which can be recognized by the apparent presence 
of oblique shock fronts in the jet itself and the ambient 



medium. These shocks form as a natural consequence of 
the nonlinear steepening of oblique sound waves as their 
amplitude becomes large. Practically all the four quanti- 
ties displayed for a given model (i.e. pressure, rest-mass 
density, Lorentz factor and internal energy) show the same 
structure, so there is no need to discuss them separately. 
Therefore we treat the pressure distribution as a repre- 
sentative quantity. The following properties of the flow 
patterns can be noticed in Figs. I4I7I 

1. The interface between jet and ambient medium forms 
a regular sinusoidal pattern. The amplitudes are all 
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Fig. 3. Pressure distribution at the end of linear phase for models BOS (left panel) and DOS (right panel). The 
continuous line corresponds to jet mass fraction equal to 0.5 and serves to distinguish jet and ambient media. 



comparable and have values of the order of O.S — 1.0i?j. 
In all cases the departure from a regular sinusoidal 
shape is rather small. 

2. The structure of oblique shocks crosses the jet inter- 
face in such a way that the deformed interface is almost 
parallel to the shock front. One should notice however 
that the oblique shock moves in the flow direction in 
the external medium and in the counter-flow direction 
in the jet medium, as can be deduced from the distri- 
bution of hot post-shock and cold pre-shock material. 
These two oppositely moving shock waves fit together 
because of the jet background flow. 

3. The highest gas pressure is always located on the 
jet symmetry plane (the brightest spots on pressure 
plots), where the crossing shocks form the familiar x- 
pattern, but the pressure enhancement on the jet inter- 
face is almost as strong as on the jet symmetry plane. 

4. The shocks are much stronger (i.e., the jumps between 
pre- and post-shock pressure larger) in the cold cases. 
In the case of hot models the jumps are smoother. This 
is a straightforward consequence of the larger ratios of 
the perturbated velocity to the sound speed in the cold 
cases. 

5. There is no significant variation of the structure of non- 
linear perturbations for different values of the Lorcntz 
factor. 

6. Model C20 displays a structure of internal waves, 
which is absent in other cases. This fine structure is 
probably connected to the distinct evolution of this 
model in the postlinear phase (see Sect. 13. H and Fig. [21). 

Therefore the differences between models just before 
the end of the saturation phase can be considered as mi- 
nor. 



4. Summary 

This is the first of a series of papers devoted to the ef- 
fects of relativistic dynamics and thermodynamics in the 
development of KH instabilities in planar relativistic jets, 
accross both linear and fully nonlinear regimes. To this 
aim, we have performed a linear stability analysis and 
high-resolution numerical simulations for the most unsta- 
ble first reflection modes in the temporal approach, for 



three different values of the jet Lorentz factor 7 (5, 10 
and 20) and a few different values of specific internal en- 
ergy of the jet matter (from 0.08 to 60. Oc^). In this paper 
we concentrate in the early stages of evolution of the KH 
instability, namely the linear and saturation phases. The 
present paper is also intended to set the theoretical and 
numerical background for the whole series of papers. 

Our simulations describe the linear regime of evolution 
of the excited eigenmodes of the different models with 
high accuracy. The growth rates of the perturbed modes 
in the vortex sheet approximation were determined with 
an average relative error of 20%. 

In all the examined cases the longitudinal velocity per- 
turbation is the first quantity that departs from linear 
growth when it reaches a value close to the speed of light 
in the jet reference frame. The reason for this saturation, 
specific to relativistic dynamics, is not so obvious in the 
reference frame of the external medium where it saturates 
at a smaller value (~ c/7, where 7 is the bulk flow Lorentz 
factor in the jet). 

The saturation phase extends from the end of the lin- 
ear phase up to the saturation of the transversal velocity 
perturbation (at approximately O.Sc in the jet reference 
frame). The saturation times for the different numerical 
models have been explained from elementary consider- 
ations, i.e. from properties of linear modes provided by 
the linear stability analysis and from the limitation of the 
transversal perturbation velocity. 

The limitation of the components of the velocity per- 
turbation at the end of the linear and saturation phases 
allows us to conclude that the relativistic nature of the 
flow appears to be responsible for the departure of the 
system from linear evolution. This behaviour is consistent 
with the one deduced by Hanasz (1995, 1997) with the aid 
of analytical methods. 

At the saturation time the perturbation structure is 
close to the structure of the initial perturbation (the one 
corresponding to the most unstable first reflection mode), 
except that the oblique sound waves forming the pertur- 
bation became steep due to their large amplitude. It is 
interesting to note that the oblique shocks are stronger 
(i.e., the pressure jumps are larger) in the colder cases. 
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Fig. 4. Snapshot around saturation of logarithmic maps of pressure, rest-mass density and specific internal energy and 
non- logarithmic Lorentz factor for model AOS. 



Our simulations, performed for the most unstable first 
reflection modes, confirm the general trends resulting from 
the linear stability analysis: the faster (larger Lorentz fac- 
tor) and colder jets have smaller growth rates. As we 
mentioned in the Introduction, Hardee et al. (1998) and 
Rosen et al. (1999) note an exception which occurs for 
the hottest jets. These jets appear to be the most sta- 
ble in their simulations (see also the simulations in Marti 
et al. 1997). They suggest that this behaviour is caused 
by the lack of appropriate perturbations to couple to the 
unstable modes. This could be partially true as fast, hot 
jets do not generate overpressured coocons that let the 
jet run directly into the nonlinear regime. However, from 
the point of view of our results, the high stability of hot 
jets may have been caused by the lack of radial resolution 
that leads to a damping in the perturbation growth rates. 
We note that the agreement between the linear stability 
analysis and numerical simulations of KH instability in 
the linear range has been achieved for a very high radial 
resolution of 400 zones /Rj, which appears to be especially 
relevant for hot jets. Finally, one should keep in mind that 
the simulations performed in the aforementioned papers 
only covered about one hundred time units, well inside 
the linear regime of the corresponding models for small 
perturbations. 

The high accuracy of our simulations in describing the 
early stages of evolution of the KH instability (as derived 
from the agreement between the computed and expected 
linear growth rates and the consistency of the saturation 
times) establishes a solid basis to study the fully nonlinear 
regime, to be done in the second paper of this series. In 
this paper (Paper II) we will show that the similarities 
found in the evolution of all the models accross the linear 
and saturation phases is lost and very different nonlinear 
evolution is found depending on the initial jet parameters. 
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Appendix A: Influence of the numerical resolution 
in the description of the linear regime of KH 
modes of relativistic flows 

Growth of the instability depends critically on the numer- 
ical viscosity of the algorithm. Hence our first aim was to 
look for suitable numerical resolutions by comparing nu- 
merical and analytical results for the linear regime. We 
performed a number of simulations based on model AOS 
changing longitudinal and radial resolution, and also the 
exponent m for the shear layer steepness (see Ean. l33| and 
I84|l . Figure IXtI shows results for several of those simula- 
tions. 

A shear layer was included in order to avoid non- 
steadiness in initial conditions given by discontinuous 
separation between the beam and the external medium. 
Therefore, if for some transversal resolution we introduced 
a shear layer which was too thin, we found a similar non- 
steady behavior. However, in order to reproduce the linear 
regime, we needed a steep enough shear layer, due to the 
fact that theory was developed for a discontinuous sepa- 
ration between both media. We can see in figure rO that 
TO = 10 models are considerably damped with respect 
to the theoretical growth, independent of the transver- 
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sal resolution. Hence resolution perpendicular to the flow 
appeared to be essential, requiring very high resolutions 
(400 zones/i?j) and thin shear layers (m = 40) with 40 to 
45 zones, i.e., an equilibrium between steepness and num- 
ber of cells. Lower transversal resolutions and/or thicker 
shear layers led to non-satisfactory results, with a slow or 
damped growth. A very low longitudinal resolution results 
in damping of growth rate, too, as can be seen from com- 
parison of m = 20 models, but Tz = 16 seems to give a 
reasonable growth rate compared to theory. This (small) 
resolution of 16 zones/ Rj along the jet was taken as a com- 
promise between accuracy and computational efficiency. 
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Fig. A.l. Linear growth of the amplitude of the pressure perturbation (in logarithmic scale) versus time (in units of 
Rj/c) for different resolutions, stands for transversal resolution, for longitudinal resolution, and m is the value 
which gives the shear layer steepness. 



